/*
 Copyright (c) 2013 MeherTJ G. All rights reserved.
 License: LGPL for personnal study or free software.
 */

#pragma once

#include <memory>

#include "zen_vector.h"
#include "zen_quaternion.h"

namespace Zen {
	struct Mat2 {
		float c[2][2];
		
		float & operator[](size_t i) {
			return (*c)[i];
		}
		float const & operator[](size_t i) const {
			return (*c)[i];
		}
	};

	struct Mat3 {
		float c[3][3];
		
		float & operator[](size_t i) {
			return (*c)[i];
		}
		float const & operator[](size_t i) const {
			return (*c)[i];
		}
	};

	struct Mat4 {
		float c[4][4];
		/*
		 // m30, m31, and m32 correspond to the translation values tx, ty, tz, respectively.
		 float m00, float m01, float m02, float m03,
		 float m10, float m11, float m12, float m13,
		 float m20, float m21, float m22, float m23,
		 float m30, float m31, float m32, float m33,
		 */
		
		float & operator[](size_t i) {
			return (*c)[i];
		}
		float const & operator[](size_t i) const {
			return (*c)[i];
		}
	};

	typedef const Mat2 cMat2;
	typedef const Mat3 cMat3;
	typedef const Mat4 cMat4;
}

inline Zen::Mat3 operator+(Zen::cMat3 & mL, Zen::cMat3 & mR);
inline Zen::Mat3 operator-(Zen::cMat3 & mL, Zen::cMat3 & mR);
inline Zen::Mat3 operator*(Zen::cMat3 & mL, Zen::cMat3 & mR);
inline Zen::Mat3 operator/(Zen::cMat3 & mL, Zen::cMat3 & mR);
inline Zen::Mat4 operator+(Zen::cMat4 & mL, Zen::cMat4 & mR);
inline Zen::Mat4 operator-(Zen::cMat4 & mL, Zen::cMat4 & mR);
inline Zen::Mat4 operator*(Zen::cMat4 & mL, Zen::cMat4 & mR);
inline Zen::Mat4 operator/(Zen::cMat4 & mL, Zen::cMat4 & mR);
inline Zen::Mat4 operator&(Zen::cMat4 & mL, Zen::cMat4 & mR);

namespace Zen {

	inline Mat2 Mat2Make(cMat3 & matrix)
	{
return {
			matrix[0], matrix[1],
			matrix[3], matrix[4]
		};
	}
	inline Mat2 Mat2Make(cMat4 & matrix)
	{
return {
			matrix[0], matrix[1],
			matrix[4], matrix[5]
		};
	}

	inline Mat3 Mat3Make(std::initializer_list<float> values)
	{
		Mat3 m3;
		float * dest = &(**m3.c);
		for(auto i : values) *dest++ = i;
		return m3;
	}

	inline Mat3 Mat3Make(cMat4 & matrix)
	{
		return {matrix[0], matrix[1], matrix[2],
			matrix[4], matrix[5], matrix[6],
			matrix[8], matrix[9], matrix[10] };
	}

	inline Mat3 Mat3MakeIdentity(float v = 1.f)
	{
		return Mat3{v, 0, 0, 0, v, 0, 0, 0, v};
	}
	
	inline Mat3 Mat3MakeWithRows(cVec3 & row0, cVec3 & row1, cVec3 & row2)
	{
		return {
			row0[0], row1[0], row2[0],
			row0[1], row1[1], row2[1],
			row0[2], row1[2], row2[2]
		};
	}
	
	inline Mat3 Mat3MakeWithColumns(cVec3 & col0, cVec3 & col1, cVec3 & col2)
	{
		return {
			col0[0], col0[1], col0[2],
			col1[0], col1[1], col1[2],
			col2[0], col2[1], col2[2] };
	}
	inline Mat3 Mat3MakeWithQuat(Vec4 const & q)
	{
		auto qua = QuatNormalize(q);
		
		float x = qua[0];
		float y = qua[1];
		float z = qua[2];
		float w = qua[3];
		
		float _2x = x + x;
		float _2y = y + y;
		float _2z = z + z;
		float _2w = w + w;
		
		return {1.0f - _2y * y - _2z * z,
			_2x * y + _2w * z,
			_2x * z - _2w * y,
			
			_2x * y - _2w * z,
			1.0f - _2x * x - _2z * z,
			_2y * z + _2w * x,
			
			_2x * z + _2w * y,
			_2y * z - _2w * x,
			1.0f - _2x * x - _2y * y };
		
	}
	
	inline Mat3 Mat3MakeScale(cVec3 & s)
	{
		return Mat3{
			s.x, 0, 0,
			0, s.y, 0,
			0, 0, s.z,
		};
	}
	
	inline Mat3 Mat3MakeRotation(float radians, cVec3 & axis)
	{
		Vec3 v = Vec3Normalize(axis);
		float cos = cosf(radians);
		float cosp = 1.0f - cos;
		float sin = sinf(radians);
		
		return {
			cos + cosp * v[0] * v[0],
			cosp * v[0] * v[1] + v[2] * sin,
			cosp * v[0] * v[2] - v[1] * sin,
			
			cosp * v[0] * v[1] - v[2] * sin,
			cos + cosp * v[1] * v[1],
			cosp * v[1] * v[2] + v[0] * sin,
			
			cosp * v[0] * v[2] + v[1] * sin,
			cosp * v[1] * v[2] - v[0] * sin,
			cos + cosp * v[2] * v[2] };
		
	}
	
	inline Mat3 Mat3MakeXRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
return { 1.0f, 0.0f, 0.0f,
			0.0f, cos, sin,
			0.0f, -sin, cos };
		
	}
	
	inline Mat3 Mat3MakeYRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
		return {cos, 0.0f, -sin,
			0.0f, 1.0f, 0.0f,
			sin, 0.0f, cos };
		
	}
	
	inline Mat3 Mat3MakeZRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
		return {
			cos, sin, 0.0f,
			-sin, cos, 0.0f,
			0.0f, 0.0f, 1.0f
		};
		
	}
	
	inline Mat2 Mat3GetMat2(cMat3 & matrix)
	{
		Mat2 m = { matrix[0], matrix[1],
			matrix[3], matrix[4] };
		return m;
	}
	
	inline Vec3 Mat3GetRow(cMat3 & matrix, int row)
	{
		Vec3 v = { matrix[row], matrix[3 + row], matrix[6 + row] };
		return v;
	}
	
	inline Vec3 Mat3GetColumn(cMat3 & matrix, int column)
	{
		Vec3 v = { matrix[column * 3 + 0], matrix[column * 3 + 1], matrix[column * 3 + 2] };
		return v;
	}
	
	inline void Mat3SetRow(Mat3 & matrix, int row, Vec3 vector)
	{
		matrix[row] = vector[0];
		matrix[row + 3] = vector[1];
		matrix[row + 6] = vector[2];
	}
	
	inline void Mat3SetColumn(Mat3 & matrix, int column, Vec3 vector)
	{
		matrix[column * 3 + 0] = vector[0];
		matrix[column * 3 + 1] = vector[1];
		matrix[column * 3 + 2] = vector[2];
	}
	
	inline Mat3 Mat3Transpose(cMat3 & o)
	{
		return {
			o[0], o[3], o[6],
			o[1], o[4], o[7],
			o[2], o[5], o[8]
		};
	}
	
	inline Mat3 Mat3Multiply(cMat3 & mL, cMat3 & mR)
	{
		Mat3 m;
		
		m[0] = mL[0] * mR[0] + mL[3] * mR[1] + mL[6] * mR[2];
		m[3] = mL[0] * mR[3] + mL[3] * mR[4] + mL[6] * mR[5];
		m[6] = mL[0] * mR[6] + mL[3] * mR[7] + mL[6] * mR[8];
		
		m[1] = mL[1] * mR[0] + mL[4] * mR[1] + mL[7] * mR[2];
		m[4] = mL[1] * mR[3] + mL[4] * mR[4] + mL[7] * mR[5];
		m[7] = mL[1] * mR[6] + mL[4] * mR[7] + mL[7] * mR[8];
		
		m[2] = mL[2] * mR[0] + mL[5] * mR[1] + mL[8] * mR[2];
		m[5] = mL[2] * mR[3] + mL[5] * mR[4] + mL[8] * mR[5];
		m[8] = mL[2] * mR[6] + mL[5] * mR[7] + mL[8] * mR[8];
		
		return m;
	}
	
	inline Mat3 Mat3Scale(cMat3 & matrix, float x, float y, float z)
	{
		return {matrix[0] * x, matrix[1] * x, matrix[2] * x,
			matrix[3] * y, matrix[4] * y, matrix[5] * y,
			matrix[6] * z, matrix[7] * z, matrix[8] * z };
	}
	
	inline Mat3 Mat3Scale(cMat3 & matrix, Vec3 scaleVector)
	{
		return {matrix[0] * scaleVector[0], matrix[1] * scaleVector[0], matrix[2] * scaleVector[0],
			matrix[3] * scaleVector[1], matrix[4] * scaleVector[1], matrix[5] * scaleVector[1],
			matrix[6] * scaleVector[2], matrix[7] * scaleVector[2], matrix[8] * scaleVector[2] };
	}
	
	inline Mat3 Mat3Scale(cMat3 & matrix, Vec4 scaleVector)
	{
		return {matrix[0] * scaleVector[0], matrix[1] * scaleVector[0], matrix[2] * scaleVector[0],
			matrix[3] * scaleVector[1], matrix[4] * scaleVector[1], matrix[5] * scaleVector[1],
			matrix[6] * scaleVector[2], matrix[7] * scaleVector[2], matrix[8] * scaleVector[2] };
	}
	
	inline Mat3 Mat3Rotate(cMat3 & matrix, float radians, cVec3 & axis)
	{
		Mat3 rm = Mat3MakeRotation(radians, axis);
		return Mat3Multiply(matrix, rm);
	}

	inline Mat3 Mat3RotateX(cMat3 & matrix, float radians)
	{
		Mat3 rm = Mat3MakeXRotation(radians);
		return Mat3Multiply(matrix, rm);
	}
	
	inline Mat3 Mat3RotateY(cMat3 & matrix, float radians)
	{
		Mat3 rm = Mat3MakeYRotation(radians);
		return Mat3Multiply(matrix, rm);
	}
	
	inline Mat3 Mat3RotateZ(cMat3 & matrix, float radians)
	{
		Mat3 rm = Mat3MakeZRotation(radians);
		return Mat3Multiply(matrix, rm);
	}
	
	inline Vec3 Mat3Multiply(cMat3 & mL, Vec3 vR)
	{
		return { mL[0] * vR[0] + mL[3] * vR[1] + mL[6] * vR[2],
			mL[1] * vR[0] + mL[4] * vR[1] + mL[7] * vR[2],
			mL[2] * vR[0] + mL[5] * vR[1] + mL[8] * vR[2] };
	}

	inline Mat4 Mat4Make(std::initializer_list<float> const & values)
	{
		Mat4 m;
		auto dest = &(**m.c);
		for(auto i : values) *dest++ = i;
		return m;
	}
	inline Mat4 Mat4MakeIdentity(float v = 1.f)
	{
		return Mat4{ v, 0, 0, 0, 0, v, 0, 0, 0, 0, v, 0, 0, 0, 0, v };
	}
	
	inline Mat4 Mat4MakeWithRows(cVec4 & row0,
								 cVec4 & row1,
								 cVec4 & row2,
								 cVec4 & row3)
	{
		return { row0[0], row1[0], row2[0], row3[0],
			row0[1], row1[1], row2[1], row3[1],
			row0[2], row1[2], row2[2], row3[2],
			row0[3], row1[3], row2[3], row3[3] };
	}
	
	inline Mat4 Mat4MakeWithColumns(cVec4 & col0,
									cVec4 & col1,
									cVec4 & col2,
									cVec4 & col3)
	{
		return {
			col0[0], col0[1], col0[2], col0[3],
			col1[0], col1[1], col1[2], col1[3],
			col2[0], col2[1], col2[2], col2[3],
			col3[0], col3[1], col3[2], col3[3] };
	}
	
	inline Mat4 Mat4MakeWithQuat(Vec4 qua)
	{
		qua = QuatNormalize(qua);
		
		float x = qua[0];
		float y = qua[1];
		float z = qua[2];
		float w = qua[3];
		
		float _2x = x + x;
		float _2y = y + y;
		float _2z = z + z;
		float _2w = w + w;
		
		return {1.0f - _2y * y - _2z * z,
			_2x * y + _2w * z,
			_2x * z - _2w * y,
			0.0f,
			_2x * y - _2w * z,
			1.0f - _2x * x - _2z * z,
			_2y * z + _2w * x,
			0.0f,
			_2x * z + _2w * y,
			_2y * z - _2w * x,
			1.0f - _2x * x - _2y * y,
			0.0f,
			0.0f,
			0.0f,
			0.0f,
			1.0f };
	}
	
	inline Mat4 Mat4MakeTranslation(cVec3 & t)
	{
		cMat4 & m = {
			1, 0, 0, 0,
			0, 1, 0, 0,
			0, 0, 1, 0,
			t.x, t.y, t.z, 1,
		};
		return m;
	}
	
	inline Mat4 Mat4MakeScale( cVec3 & s)
	{
		return {
			s.x, 0, 0, 0,
			0, s.y, 0, 0,
			0, 0, s.z, 0,
			0, 0, 0, 1,
		};
	}

	/// Rotation begin
	inline Mat4 Mat4MakeRotation(float radians, cVec3 & r)
	{
		Vec3 axis = Vec3Normalize(r);
		float ct = cosf(radians);
		float st = sinf(radians);
		float ci = 1 - ct;
		float x = axis.x;
		float y = axis.y;
		float z = axis.z;
		
		return Mat4{
			ct + x * x * ci, y * x * ci + z * st, z * x * ci - y * st, 0,
			x * y * ci - z * st, ct + y * y * ci, z * y * ci + x * st, 0,
			x * z * ci + y * st, y * z * ci - x * st, ct + z * z * ci, 0,
			0, 0, 0, 1,
		};
	}
	
	inline Mat4 Mat4MakeXRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
		return {
			1.0f, 0.0f, 0.0f, 0.0f,
			0.0f, cos, sin, 0.0f,
			0.0f, -sin, cos, 0.0f,
			0.0f, 0.0f, 0.0f, 1.0f
		};
	}
	
	inline Mat4 Mat4MakeYRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
		return {
			cos, 0.0f, -sin, 0.0f,
			0.0f, 1.0f, 0.0f, 0.0f,
			sin, 0.0f, cos, 0.0f,
			0.0f, 0.0f, 0.0f, 1.0f
		};
		
	}
	
	inline Mat4 Mat4MakeZRotation(float radians)
	{
		float cos = cosf(radians);
		float sin = sinf(radians);
		
		return {
			cos, sin, 0.0f, 0.0f,
			-sin, cos, 0.0f, 0.0f,
			0.0f, 0.0f, 1.0f, 0.0f,
			0.0f, 0.0f, 0.0f, 1.0f
		};
	}
	
	inline Mat4 Mat4MakePerspective(float fovyRadians, float aspect, float /*nearZ*/ n, float /*farZ*/ f)
	{
		float ys = 1 / tanf(fovyRadians * 0.5);
		float xs = ys / aspect;
		float dz = n - f;
		float zs = (n+f) / dz;
		
		return {
			xs, 0, 0, 0,
			0, ys, 0, 0,
			0, 0, zs, -1,
			0, 0, 2 * f * n / dz, 0// 2 * f * n / dz - zs * n, +n,
		};
		
		/* move nearZ to zero.
		 return {
		 xs, 0, 0, 0,
		 0, ys, 0, 0,
		 0, 0, zs, -1,
		 0, 0, 2 * f * n / dz + 2 * f * n / dz - zs * n, +n,
		 };
		 */
	}
	
	inline Mat4 Mat4MakeFrustum(float left, float right,
								float bottom, float top,
								float nearZ, float farZ)
	{
		float ral = right + left;
		float rsl = right - left;
		float tsb = top - bottom;
		float tab = top + bottom;
		float fan = farZ + nearZ;
		float fsn = farZ - nearZ;
		
		return {
			2.0f * nearZ / rsl, 0.0f, 0.0f, 0.0f,
			0.0f, 2.0f * nearZ / tsb, 0.0f, 0.0f,
			ral / rsl, tab / tsb, -fan / fsn, -1.0f,
			0.0f, 0.0f, (-2.0f * farZ * nearZ) / fsn, 0.0f
		};
	}
	
	inline Mat4 Mat4MakeOrtho(float left, float right,
							  float bottom, float top,
							  float nearZ, float farZ)
	{
		float ral = right + left;
		float rsl = right - left;
		float tab = top + bottom;
		float tsb = top - bottom;
		float fan = farZ + nearZ;
		float fsn = farZ - nearZ;
		
		return {
			2.0f / rsl, 0.0f, 0.0f, 0.0f,
			0.0f, 2.0f / tsb, 0.0f, 0.0f,
			0.0f, 0.0f, -2.0f / fsn, 0.0f,
			-ral / rsl, -tab / tsb, -fan / fsn, 1.0f
		};
	}
	
	inline Mat4 Mat4MakeLookAt(cVec3 & ev, cVec3 & cv, cVec3 & uv)
	{
		Vec3 n = Vec3Normalize(ev - cv);
		Vec3 u = Vec3Normalize(Vec3CrossProduct(uv, n));
		Vec3 v = Vec3CrossProduct(n, u);
		
		return {u[0], v[0], n[0], 0.0f,
			u[1], v[1], n[1], 0.0f,
			u[2], v[2], n[2], 0.0f,
			Vec3DotProduct(-(u), ev),
			Vec3DotProduct(-(v), ev),
			Vec3DotProduct(-(n), ev),
			1.0f };
		
	}

	
	inline Vec4 Mat4GetRow(cMat4 & matrix, int row)
	{
		Vec4 v = { matrix[row], matrix[4 + row], matrix[8 + row], matrix[12 + row] };
		return v;
	}
	
	inline Vec4 Mat4GetColumn(cMat4 & matrix, int column)
	{
		Vec4 v = { matrix[column * 4 + 0], matrix[column * 4 + 1], matrix[column * 4 + 2], matrix[column * 4 + 3] };
		return v;
	}
	
	inline void Mat4SetRow(Mat4 & matrix, int row, cVec4 & v)
	{
		matrix[row] = v[0];
		matrix[row+4] = v[1];
		matrix[row+8] = v[2];
		matrix[row+12] = v[3];
	}
	
	inline void Mat4SetColumn(Mat4 & matrix, int column, cVec4 & v)
	{
		matrix[column*4+0] = v[0];
		matrix[column*4+1] = v[1];
		matrix[column*4+2] = v[2];
		matrix[column*4+3] = v[3];
	}

	/// multiply
	inline Vec4 Mat4Multiply(cMat4 & mL, cVec4 & vR)
	{
		return {
			mL[0] * vR[0] + mL[4] * vR[1] + mL[8] * vR[2] + mL[12] * vR[3],
			mL[1] * vR[0] + mL[5] * vR[1] + mL[9] * vR[2] + mL[13] * vR[3],
			mL[2] * vR[0] + mL[6] * vR[1] + mL[10] * vR[2] + mL[14] * vR[3],
			mL[3] * vR[0] + mL[7] * vR[1] + mL[11] * vR[2] + mL[15] * vR[3]
		};
	}

	inline Vec3 Mat4Multiply(cMat4 & mL, cVec3 & vR)
	{
		Vec4 v4 = Mat4Multiply(mL, Vec4Make(vR, 0.f));
		return Vec3Make(v4[0], v4[1], v4[2]);
	}

	inline Vec3 Mat4MultiplyAndProject(cMat4 & mL, cVec3 & vR)
	{
		Vec4 v4 = Mat4Multiply(mL, Vec4Make(vR, 1.0f));
		return Vec3Make(v4[0], v4[1], v4[2]) * (1.0f / v4[3]);
	}

	inline Mat4 Mat4Multiply(cMat4 & mL, cMat4 & mR)
	{
		Mat4 m;

		m[0] = mL[0] * mR[0] + mL[4] * mR[1] + mL[8] * mR[2] + mL[12] * mR[3];
		m[4] = mL[0] * mR[4] + mL[4] * mR[5] + mL[8] * mR[6] + mL[12] * mR[7];
		m[8] = mL[0] * mR[8] + mL[4] * mR[9] + mL[8] * mR[10] + mL[12] * mR[11];
		m[12] = mL[0] * mR[12] + mL[4] * mR[13] + mL[8] * mR[14] + mL[12] * mR[15];

		m[1] = mL[1] * mR[0] + mL[5] * mR[1] + mL[9] * mR[2] + mL[13] * mR[3];
		m[5] = mL[1] * mR[4] + mL[5] * mR[5] + mL[9] * mR[6] + mL[13] * mR[7];
		m[9] = mL[1] * mR[8] + mL[5] * mR[9] + mL[9] * mR[10] + mL[13] * mR[11];
		m[13] = mL[1] * mR[12] + mL[5] * mR[13] + mL[9] * mR[14] + mL[13] * mR[15];

		m[2] = mL[2] * mR[0] + mL[6] * mR[1] + mL[10] * mR[2] + mL[14] * mR[3];
		m[6] = mL[2] * mR[4] + mL[6] * mR[5] + mL[10] * mR[6] + mL[14] * mR[7];
		m[10] = mL[2] * mR[8] + mL[6] * mR[9] + mL[10] * mR[10] + mL[14] * mR[11];
		m[14] = mL[2] * mR[12] + mL[6] * mR[13] + mL[10] * mR[14] + mL[14] * mR[15];

		m[3] = mL[3] * mR[0] + mL[7] * mR[1] + mL[11] * mR[2] + mL[15] * mR[3];
		m[7] = mL[3] * mR[4] + mL[7] * mR[5] + mL[11] * mR[6] + mL[15] * mR[7];
		m[11] = mL[3] * mR[8] + mL[7] * mR[9] + mL[11] * mR[10] + mL[15] * mR[11];
		m[15] = mL[3] * mR[12] + mL[7] * mR[13] + mL[11] * mR[14] + mL[15] * mR[15];

		return m;
	}
	/// multiply end

	
	inline Mat4 Mat4Transpose(cMat4 & m)
	{
		return {
			m[0], m[4], m[8], m[12],
			m[1], m[5], m[9], m[13],
			m[2], m[6], m[10], m[14],
			m[3], m[7], m[11], m[15]
		};
	}

	inline Mat4 Mat4Translate(cMat4 & matrix, cVec3 & t)
	{
		return {
			matrix[0], matrix[1], matrix[2], matrix[3],
			matrix[4], matrix[5], matrix[6], matrix[7],
			matrix[8], matrix[9], matrix[10], matrix[11],
			matrix[0] * t[0] + matrix[4] * t[1] + matrix[8] * t[2] + matrix[12],
			matrix[1] * t[0] + matrix[5] * t[1] + matrix[9] * t[2] + matrix[13],
			matrix[2] * t[0] + matrix[6] * t[1] + matrix[10] * t[2] + matrix[14],
			matrix[3] * t[0] + matrix[7] * t[1] + matrix[11] * t[2] + matrix[15]
		};
	}
	
	inline Mat4 Mat4Translate(cMat4 & matrix, cVec4 & t)
	{
		return {
			matrix[0], matrix[1], matrix[2], matrix[3],
			matrix[4], matrix[5], matrix[6], matrix[7],
			matrix[8], matrix[9], matrix[10], matrix[11],
			matrix[0] * t[0] + matrix[4] * t[1] + matrix[8] * t[2] + matrix[12],
			matrix[1] * t[0] + matrix[5] * t[1] + matrix[9] * t[2] + matrix[13],
			matrix[2] * t[0] + matrix[6] * t[1] + matrix[10] * t[2] + matrix[14],
			matrix[3] * t[0] + matrix[7] * t[1] + matrix[11] * t[2] + matrix[15]
		};
	}

	inline Mat4 Mat4Scale(cMat4 & matrix, cVec3 & scale)
	{
		return {
			matrix[ 0] * scale[0],
			matrix[ 1] * scale[0],
			matrix[ 2] * scale[0],
			matrix[ 3] * scale[0],
			matrix[ 4] * scale[1],
			matrix[ 5] * scale[1],
			matrix[ 6] * scale[1],
			matrix[ 7] * scale[1],
			matrix[ 8] * scale[2],
			matrix[ 9] * scale[2],
			matrix[10] * scale[2],
			matrix[11] * scale[2],
			matrix[12], matrix[13], matrix[14], matrix[15]
		};
	}
	
	inline Mat4 Mat4Rotate(cMat4 & matrix, float radians, Vec3 axis)
	{
		Mat4 rm = Mat4MakeRotation(radians, axis);
		return Mat4Multiply(matrix, rm);
	}
	
	inline Mat4 Mat4RotateX(cMat4 & matrix, float radians)
	{
		Mat4 rm = Mat4MakeXRotation(radians);
		return Mat4Multiply(matrix, rm);
	}
	
	inline Mat4 Mat4RotateY(cMat4 & matrix, float radians)
	{
		Mat4 rm = Mat4MakeYRotation(radians);
		return Mat4Multiply(matrix, rm);
	}
	
	inline Mat4 Mat4RotateZ(cMat4 & matrix, float radians)
	{
		Mat4 rm = Mat4MakeZRotation(radians);
		return Mat4Multiply(matrix, rm);
	}
}

inline Zen::Mat3 operator+(Zen::cMat3 & mL, Zen::cMat3 & mR)
{
	return {
		mL[0] + mR[0],
		mL[1] + mR[1],
		mL[2] + mR[2],
		mL[3] + mR[3],
		mL[4] + mR[4],
		mL[5] + mR[5],
		mL[6] + mR[6],
		mL[7] + mR[7],
		mL[8] + mR[8],
	};
}

inline Zen::Mat3 operator-(Zen::cMat3 & mL, Zen::cMat3 & mR)
{
	return {
		mL[0] - mR[0],
		mL[1] - mR[1],
		mL[2] - mR[2],
		mL[3] - mR[3],
		mL[4] - mR[4],
		mL[5] - mR[5],
		mL[6] - mR[6],
		mL[7] - mR[7],
		mL[8] - mR[8],
	};
}

inline Zen::Mat3 operator*(Zen::cMat3 & mL, Zen::cMat3 & mR)
{
	return {
		mL[0] * mR[0],
		mL[1] * mR[1],
		mL[2] * mR[2],
		mL[3] * mR[3],
		mL[4] * mR[4],
		mL[5] * mR[5],
		mL[6] * mR[6],
		mL[7] * mR[7],
		mL[8] * mR[8],
	};
}

inline Zen::Mat3 operator/(Zen::cMat3 & mL, Zen::cMat3 & mR)
{
	return {
		mL[0] / mR[0],
		mL[1] / mR[1],
		mL[2] / mR[2],
		mL[3] / mR[3],
		mL[4] / mR[4],
		mL[5] / mR[5],
		mL[6] / mR[6],
		mL[7] / mR[7],
		mL[8] / mR[8],
	};
}

inline Zen::Mat4 operator+(Zen::cMat4 & mL, Zen::cMat4 & mR)
{
	return {
		mL[ 0] + mR[ 0],
		mL[ 1] + mR[ 1],
		mL[ 2] + mR[ 2],
		mL[ 3] + mR[ 3],
		mL[ 4] + mR[ 4],
		mL[ 5] + mR[ 5],
		mL[ 6] + mR[ 6],
		mL[ 7] + mR[ 7],
		mL[ 8] + mR[ 8],
		mL[ 9] + mR[ 9],
		mL[10] + mR[10],
		mL[11] + mR[11],
		mL[12] + mR[12],
		mL[13] + mR[13],
		mL[14] + mR[14],
		mL[15] + mR[15],
	};
}

inline Zen::Mat4 operator-(Zen::cMat4 & mL, Zen::cMat4 & mR)
{
	return {
		mL[ 0] - mR[ 0],
		mL[ 1] - mR[ 1],
		mL[ 2] - mR[ 2],
		mL[ 3] - mR[ 3],
		mL[ 4] - mR[ 4],
		mL[ 5] - mR[ 5],
		mL[ 6] - mR[ 6],
		mL[ 7] - mR[ 7],
		mL[ 8] - mR[ 8],
		mL[ 9] - mR[ 9],
		mL[10] - mR[10],
		mL[11] - mR[11],
		mL[12] - mR[12],
		mL[13] - mR[13],
		mL[14] - mR[14],
		mL[15] - mR[15],
	};
}

inline Zen::Mat4 operator*(Zen::cMat4 & mL, Zen::cMat4 & mR)
{
	return {
		mL[ 0] * mR[ 0],
		mL[ 1] * mR[ 1],
		mL[ 2] * mR[ 2],
		mL[ 3] * mR[ 3],

		mL[ 4] * mR[ 4],
		mL[ 5] * mR[ 5],
		mL[ 6] * mR[ 6],
		mL[ 7] * mR[ 7],

		mL[ 8] * mR[ 8],
		mL[ 9] * mR[ 9],
		mL[10] * mR[10],
		mL[11] * mR[11],

		mL[12] * mR[12],
		mL[13] * mR[13],
		mL[14] * mR[14],
		mL[15] * mR[15],
	};
}

inline Zen::Mat4 operator/(Zen::cMat4 & mL, Zen::cMat4 & mR)
{
	return {
		mL[ 0] / mR[ 0],
		mL[ 1] / mR[ 1],
		mL[ 2] / mR[ 2],
		mL[ 3] / mR[ 3],
		mL[ 4] / mR[ 4],
		mL[ 5] / mR[ 5],
		mL[ 6] / mR[ 6],
		mL[ 7] / mR[ 7],
		mL[ 8] / mR[ 8],
		mL[ 9] / mR[ 9],
		mL[10] / mR[10],
		mL[11] / mR[11],
		mL[12] / mR[12],
		mL[13] / mR[13],
		mL[14] / mR[14],
		mL[15] / mR[15],
	};
}

inline Zen::Mat4 operator &(Zen::cMat4 & mL, Zen::cMat4 & mR)
{
	return Zen::Mat4Multiply(mL, mR);
}

inline Zen::Vec4 operator &(Zen::cMat4 & mL, Zen::cVec4 & vR)
{
	return Zen::Mat4Multiply(mL, vR);
}

inline std::ostream & operator << (std::ostream & o, Zen::cMat3 & m)
{
	for(int i = 0; i < 3; ++i)
	{
		o << m.c[i][0];
		for(int j = 1; j < 3; ++j)
		{
			o << ',' << m.c[i][j];
		}
		o << '\n';
	}
	return o;
}

inline std::ostream & operator << (std::ostream & o, Zen::cMat4 & m)
{
	using namespace std;
	for(int i = 0; i < 4; ++i)
	{
		o << std::to_string(m.c[i][0]);
		for(int j = 1; j < 4; ++j)
		{
			o << ',' << std::to_string(m.c[i][j]);
		}
		o << '\n';
	}
	return o;
}
